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BACKGROUND OF THE INVENTION 

1 . Field of the Invention 

The field of the invention is seismic prospecting. The invention relates to seismic 
imagery and more particularly to estimating the illumination fold in the migrated 3D domain 
necessary to obtain a high resolution image of the underground structure at true amplitude. 

2. Description of Related Art 

Seismic prospecting generally consists of emitting seismic waves into the sub-soil using 
one or several seismic sources, and recording seismic data corresponding to seismic waves 
reflected on geological interfaces in the sub-soil (also called reflectors) on the ground surface, as 
a function of time, using receivers (also called geophones or hydrophones depending on whether 
the survey is being made on land or at sea) and then processing these data to extract useful 
information about the geology of the sub-soil. 
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Data processing (or seismic traces) recorded by receivers are useful particularly for 
designing images with information about the geology of the sub-soil. 

BRffiF SUMMARY OF THE INVENTION 
5 The objective of high resolution seismic imagery with true amplitude is to provide an 

image of the sub-soil that is not disturbed, particularly by wave field propagation effects in the 
sub-soil and by acquisition effects (we then talk about acquisition patterns). 

Acoustic reflectivity of the sub-soil is represented as an image using a so-called 
migration method to restore correct forms of geological interfaces. To achieve this, seismic data 
10 are converted from the domain in which they are acquired to the migrated domain in which 
geological interfaces are represented in their genuine locations. 

In the migrated domain, the vertical scale may be time (we then talk about time 
migration) or depth (we then talk about depth migration). 

Migration is usually done before summation (we then talk about pre-stack migration), 
15 with direct mapping of the seismic signal from seismic traces towards the migrated domain. 

The precise identification of the sub-soil geology at a particular point is directly related 
to the capability of a seismic emission to illuminate this point. If a point in the sub-soil is highly 
illuminated, in all likelihood the seismic survey can evaluate the geology of the sub-soil at this 
point. On the other hand, if a point of the sub-soil is not well illuminated, it is quite possible that 
20 the seismic survey will provide an imprecise or even incorrect evaluation of the sub-soil geology 
at this point. 

In order to obtain a precise representation of the sub-soil, the true amplitudes of the 
reflection of the seismic wave on a reflector must be represented perfectly. 

Therefore, it is necessary to know precisely if a change in the amplitude at a specific 
25 point in the sub-soil is the result of a change in the geology of the sub-soil at this point or if it is 
the result of a particular illumination at this point. 

But reflectors are usually illuminated non-uniformly. In particular, the illumination 
depends on: 

- the configuration of seismic data acquisition; 
30 - lateral and vertical variations of the propagation velocity of seismic waves in the 

sub-soil. 
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In order to obtain a high resolution image of the reflectivity of the sub-soil, it is necessary 
to make an amplitude correction to take account of non-uniformity of illumination of reflectors 
when migrating seismic reflections. 

As has already been described above, this correction (done using a weighting factor) is 
5 necessary because: 

- the reflectivity of the sub-soil cannot be precisely imaged without correcting 
amplitude distortions related to propagation of the wave field in the sub-soil (spherical 
divergence, absorption, directivity, etc.); 

- to improve the resolution of the reflectivity image, it is important to take account 
10 of the variability of the migration illumination fold for each inclination of a reflector used during 

the migration. 

For a given image point in the migrated domain, the illumination fold corresponds to the 
number of source-receiver pairs for which contributions for building up the amplitude of the 
reflectivity image interfere constructively (stationary phase condition). 
15 The document "Bleistein, N. 1987, On the Imaging of Reflectors in the Earth: 

Geophysics Soc. of Expl. Geophysics, 931-942" proposes a method used to obtain an amplitude 
weighting factor for a regular acquisition configuration, capable of making a high resolution 
representation of the reflectivity of the sub-soil in true amplitude. 

This is an "inverse" method composed of a modelling step (modelling of seismic traces 
20 from a reflectivity model) followed by a migration step (from modelled seismic traces to imaged 
reflectivity), the weighting factor then being chosen so as to obtain perfect agreement between 
the imaged reflectivity and the reflectivity model. 

By following this "inverse" method proposed by Bleistein and as shown in Figure 1, it 
can be seen that: 

25 - for a given acquisition configuration (distribution of (source S, receiver R) pairs), 

- for a given image point x in the migrated domain, 

- for a given reflector inclination (indicated by a dip vector p ), 

the illumination fold | (jc, p) is determined by the number of source-receiver pairs such 
that their migration contributions constructively interfere with the image of the dip reflector p at 
30 image point x . 
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Expressions developed by Bleistein for the illumination fold | (x, p) include the effect 
of the entire geometric configuration of the acquisition in a single step. 

Consequently, the method proposed by Bleistein can only be used for a regular and 
parametric distribution of all source-receiver pairs. 
5 Consequently, this restrictive application condition limits its use to a small number of 

acquisition configurations. 

Therefore for acquisition configurations with irregular geometries, it is necessary to 
define individual illumination folds | (3c, p, s, r), for each (source s , receiver r ) pair involved in 
the migration. Obviously, the said individual ratios contribute to the total illumination fold 
10 I (x, p) according to I (I, p) = 

^ / (3c,p,J,r). 

sr 

However, the study carried out by Bleistein is based on a number of developments that 
were intended to estimate the illumination fold in the migrated domain for any geometric 
acquisition configuration, and for any distribution of the source receiver pairs. 
15 For years, various approximations have been proposed for the seismic wave that is the 

key element for estimating the stationary phase condition (constructive interference condition). 

The first developments consisted of modeling specular reflection points along a horizon 
(2D map) selected in the migrated domain (3D) by ray tracing for each of the (source J, 
receiver r) pairs. 

20 This reflection modelling step determines the following quantities: 

location of reflection points (s, r) along the selected horizon; 
- the specular reflection travel time (x; J, r) ; 

The diffraction travel time td along the selected horizon can then be estimated by a 
development of the reflection travel time , around the specular reflection point x^ (?, r ) in the 
25 Taylor series (determination of first and second derivatives Vj^ (x^ ) and A^ J^ (x^ ) ). 

Figure 2 illustrates this procedure during which the diffraction travel time (x ) at image 
point X is estimated using the reflection travel time tr(x^)and more particularly by its 
development around the specular reflection point x^ in the Taylor series. 

30 
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Summation of all wavelet mappings obtained for each source-receiver pair then provides 
a means of estimating the illumination fold in the migrated domain along the selected horizon. 

By modelling the seismic wavelet according to a pulse peak, with infinite frequency 
band, document "Cain, G., Cambois, G., Gehin, M and Hall, R., 1998 Reducing risk in seismic 
5 acquisition and interpretation of complex targets using a Gocad-based 3-D modeling tool, 68th 
, Ann. Intemat. Mtg, Soc. of Expl. Geophys., 2072-2075" has shown that an impact map can be 
obtained, and in this case only the specular contribution on a selected horizon is actually taken 
into account. 

Considering a realistic approximation of a seismic wavelet with limited band, the 
10 following documents: 

- "Schneider, W., Jr and Winbow, G., 1999, Efficient and accurate modelling of 3D 
seismic illumination, 69th Ann. Intemat. Mtg, Soc. Of ExpL Geophys, Expanded Abstract, 633- 
636"; 

- "Laurain, R. and Vinje, V. 2001, Pre-stack migration and illumination maps, 71th 
15 Ann. Intemat., Mtg: Soc. Of Expl. Geophys, Expanded Abstract, 929-932"; have shown that it is 

possible to improve estimating the illumination fold along selected horizons. 

These approaches, in particular based on scanning of selected horizons, have been found 
to be reliable and precise. However, their efficiency is limited particularly due to the fact that: 

- compensation for the illumination fold is limited to a finite number of selected 
20 horizons and therefore, unlike the method developed by Bleistein, this method is not completely 

3D; 

- for each (source J, receiver r) pair involved in the migration, a dedicated 
modelling step must be carried out after the migration. However, this modelling step to 
determine the(jc^, V^/^,A^^f^) quantities mentioned previously by ray tracing along each of the 

25 selected horizons is particularly expensive in terms of time. 

In order to overcome these limitations, and in order to approximate a 3D approach, a two- 
step solution has been proposed to estimate and then to compensate for the illumination ratio, 
and particularly in documents: 

"Albertin, U., Bloor, R., Beasley, C, Chang, W., Jaramillo, H., Mobley, E. and 
30 Yingst, D., 1999, Aspects of true amplitude migration, 69th Ann. Intemat. Mtg: Soc. Of Expl. 
Geophys, Expanded Abstract, 1247-1250"; 
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- "Audebert, F., Froidevaux, P. Huard, I., Nicoletis, L. and Svay-Lucas, J., 2000, A 
multi angle toolbox for restored amplitude images and AVA gathers, 70th Ann. Intemat. Mtg: 
Soc. Of Expl. Geophys, Expanded Abstract, 1020-1023". 

The first step consists of making an analysis of the multi-dip illumination fold using 
5 predetermined tables for Kirchhoff migration for diffraction travel times t^(x;s,r), in the 
migrated domain. 

In the second step, the inverse values of the estimated illumination folds are used to make 
the weighted sunmiation of seismic data. 

For each image point x and for an entire reflector dip range p , it has been proposed to 
10 increment an illumination occurrence counter I(x,p) every time that a migration contribution 
reaches the image point x with the dip component p . 

In a manner similar to that proposed by Cain et al (see above), the criterion used to 
increase the illumination occurrence counter takes account of a pulse model of the seismic 
wavelet, with infinite frequency band. 
15 This high frequency approximation is such that for a given dip component, only specular 

migration contributions corresponding to the said dip component are used for estimating the 
illumination. 

And there are a number of limitations to this high frequency approximation: 

- the estimated illumination is biased since the limited band nature of the seismic 
20 signal is ignored; 

- in order to avoid division by zeros, the illumination folds are stored by dip 
classes. However, the size of each class is a critical element. If the class size is too small, risks of 
numerical instability are very high; if the class size is too large, the precision of the estimated 
illumination coverage is low. 

25 The invention is designed to provide a method for estimating the seismic illumination 

fold in the migrated 3D domain that is not subject to these limitations. 

More precisely, the purpose of the invention is to propose a method that takes account of 
the limited band nature of the seismic signal in order to estimate more realistic illumination 
folds, and which at the same time solves the problem of the class size. 

30 
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Consequently, according to a first aspect of the invention, it divulges a method for 
estimating the seismic illumination fold | (jc, p) in the migrated 3D domain at an image point x 
for a dip of vector p , characterized in that it estimates the illumination fold I (I, p, 5 , r ) for each 
(source?, receiverr ) pair in the seismic survey, by applying the following steps: 
5 - determination of the reflection travel time t^(x^(p);sr) from the source s to the 

specular reflection point 3c^on the plane reflector passing through the image point jc and 
perpendicular to the dip vector p and then return to the reflector r ; 

- starting from the diffraction travel time tj{x;J,r) from the source sto the said 
image point x and then return to the reflector r ; 

10 ' - incrementing the said illumination fold I (x, p, J, r) related to the said (source?, 
receiverr) pair as a function of the difference between the diffraction travel time /^(x; J,F) and 
the reflection travel time t^(x^ (p);Jr) . 

Preferred but non-limitative aspects of the method according to the first aspect of the 
invention are as follows: 

15 - the method also comprises a step to summate each illumination fold 

I(jc,p,J,r)related to a (source?, receiverr) pair so as to determine the total illumination fold. 

Ux,p)=Y, l(x,pJJX 

- during the increment step, the illumination fold I(jc,p,?,r)may be incremented 
using an increment function i(r^,r^;?,r) using the equation I(3c,^,?,r)= I(jc,p,?,r) + 

20 i(f^,f,;J,r), 

- the increment function i may be a function of the seismic wavelet s(t), 

- the increment function i may be expressed as a function of the derivative of the 
seismic wavelet with respect to time, 

- it is envisaged that an a priori correction w(3c, J,r) of the illumination fold could 
25 be taken into account, which usually has to be considered during migration, 

- the determination step includes the second order Taylor series development of the 
diffraction travel time (x; ?, r ) around the image point x , 
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- the specular reflection point (p) is determined along the length of the said 
reflector such that the diffraction travel time at the said specular reflection point 

x^(p) is stationary, 

5 - the determination step uses isochronic migration maps r^(jc;5',r) specified for 

each (source J, receiver r) pair involved in the migration and each image point 3c in the 
migrated 3D domain, 

- the seismic illumination fold I(3c,p) in the migrated 3D domain is estimated 

during the Kirchhoff summation migration of seismic data recorded during the 3D seismic 
10 prospecting. 

According to another aspect, the invention proposes a method for correction of seismic 
data amplitudes recorded during a 3D seismic survey in order to compensate for the effect of 
non-uniform illumination of sub-soil reflectors, characterized in that it comprises steps consisting 
of: 

15 - estimating the illumination fold l(x,p) using the method according to the first 

aspect of the invention, 

using the inverse I * (J,^)of the said ratio as a weighting factor to be applied to 
each of the said seismic data amplitudes. 

According to another aspect, the invention proposes a method for selection of acquisition 
20 geometry among a plurality of acquisition geometries as a function of the target of 3D seismic 
prospecting, characterised in that it consists of the following steps: 

- determining the illumination fold l(x,p) by the method according to the first 
aspect of the invention, for each of the acquisition geometries considered, 

- selecting the acquisition geometry providing the optimum illumination fold as a 
25 function of the target. 
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BRffiF DESCRIPTION OF THE SEVERAL VIEWS OF THE DRAWINGS 

Other aspects, purposes and advantages of the invention will become clearer after reading 
the following detailed description of a preferred embodiment of it, given as a non-limitative 
5 example with reference to the attached drawings on which: 

- Figure 1, already commented, illustrates the constructive or non-constructive 

contribution as a function of s(t^ - f J , of the seismic trace migrated to the image of the reflector 
with a dip along 3c ; 

Figure 2, already conmiented, illustrates the estimate made according to the state 
10 of the art for the diffraction travel time at the image point (x) starting from a Taylor series 

development of the reflection travel time at the specular reflection point (x^ ) ; 

- Figure 3 illustrates the method according to the invention according to which the 
specular reflection point jcis estimated, together with the reflection travel time (J J starting 
from a Taylor series development of the diffraction travel time at the image point (x) . 

15 - Figure 4 is a flowchart diagrammatically showing the main steps in the method 

according to the invention. 

DETAILED DESCRIPTION OF THE INVENTION 

The invention proposes a method of estimating the illumination fold in the migrated 
20 domain I(jc,p;J,r) for any migrated position J, regardless of its dip component p and 
regardless of the distribution of seismic traces, using only the most basic quantities involved in 
migration by a Kirchhoff summation, namely: 

- the coordinates r ) of source-receiver pairs, where s represents the position of 

a source 5 and r represents the position of a receiver R;* 
25 - migration isochronic maps td (x; J, r ) for each of the (source J , receiver r ) pairs 

involved in the migration. 

The invention is equally applicable to imagery in time and imagery in depth. 

The invention divulges a method of the type including a modelling step followed by a 
migration step. 
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Concerning the modelling step of the reflection point x and the reflection travel time 
f^(jc^) , the method according to the invention applies the stationary phase condition to a Taylor 
series development of diffraction travel times (jc) . 

The said development is a second order development. 

Conceming the migration step, the method according to the invention uses isochronic 
migration maps available for the migration conventionally done by a Kirchhoff summation. 

For each of the (source s , receiver r ) pairs, the illumination fold l{x,p\sj) is 
updated with reference to Figure 3 using the following equation: 

I (J, p\ J, r ) = I (x , p\ r ) + 1 {td , tr\ J, r) Equation ( 1 ) 

where 

- 3c is the image point in the migrated domain x -{x^,X2,x^Y \ 

- p is the dip vector representing the selected dip component, p^iPx.Pi'^Pif y 

- s defines the coordinates of the source, s = (S|, S2, s^f 

- r defines the coordinates of the receiver, r = (ri, xi, i^)^ 

- td{x\s,r) is the diffraction travel time from the source J, as far as the migrated 
position X , and retum as far as the receiver f (path indicated by on Figure 3 and shown in 
solid lines); 

~ ^r(^r(P)'^'^) is the reflection travel time from the source ?as far as point J of 
specular reflection on the plane reflector passing through image point x and perpendicular to the 
dip vector p , and retum as far as the receiver r (path indicated by r^on Figure 3 and shown in 
dashed lines); 

- l{td,tr\sj) is an increment function for which the expression will be described 
in detail below. 

As already mentioned above, the method according to the invention advantageously uses 
isochronic migration maps t^(x;s,r) involved in the migration by Kirchhoff summation, for 

each of the (source s , receiver r ) pairs involved in the migration. These maps are "input data" 
to the method according to the invention. 
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It also uses first and second derivatives V^f^(jc; J,r),A^^r^(jc; J,r)of migration 

isochrones when they are available. If not, the said derivatives are not input data, but may be 
calculated as part of the method according to the invention. 

Figure 4 is a flow chart showing the different steps in the method according to a preferred 
5 embodiment of the invention. 

The said method consists of carrying out a step 1 used for each image point x in the 
migrated 3D domain, for each of the (source s , receiver r ) pairs involved in the migration. 

The said step 1 includes two operations 2, 3, for each dip component p selected for the 
illumination calculation I (x, p\ J, r) . 
10 The first operation 2 consists of determining the reflection travel time {x^ (p); sr) to the 

specular reflection point x^(p)on the plane reflector passing through the image point x and 
perpendicular to the dip vector p . 

Operation 2 is determined from the diffraction travel time td(x;JJ) known from the 
isochronic migration maps. 
15 The second operation 3 consists of incrementing the illumination fold in the migrated 

domain I(jc,p;J,r)of an increment function I(r^,r^;?,r) the said increment function taking 
account of: 

- the difference J,r)— r^(jc^(p); Jr) } between the travel times mentioned 

above; 

20 - the seismic wavelet s(t); 

- and possibly the migration amplitude term w(x; J, r ) 

At the end of operations 2 and 3, the illumination amplitude ratio l(x,p;sj) is known 
for a given image point x and a dip component p . 

By repeating the operations mentioned above for each of the (source s , receiver r ) pairs, 
25 the final step is to determine the illumination ratio in the migrated domain I (jc, p) 

This is done by adding the relative illumination to each pair (source J , receiver r) 
according to I (x, p) = ^ /(J, p; J, r) 
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The reflection point and the reflection travel time are determined (during 
operation 2 mentioned above) starting from a stationary phase analysis, by a second order 
Taylor series development of the diffraction travel time td around the image point x . 

The following expressions for jc^and (equations (3) and (4)) are obtained from 
this analysis (see equation (2) below): 

t, (x ; r ) = (x; J, r) + ( V/^ y {x-x^) + \l2{x-x^y {V ^jXx-x^) + ... 

Equation (2) 

jc^ (p) = jc - M.F^ • b Equation (3) 

tr (^r {p)\sj)^t^{x,sj)-\l2'b\¥^'b Equation (4) 

where, in equations (3) and (4). 

- M is a (3 X 2) matrix described as two vectors extending along the length 
of the reflection plane. Consequently, these two vectors are perpendicular to the dip 
vector p , 

- Z? is a (2 X 1) vector of the second order derivatives of the diffraction 
travel time along the reflection plane, 

F = M.^ *^^^'') ^ Equation (5) 

- F is a (2 X 2) matrix of first order derivatives of the diffraction travel time 
along the reflection plane, 

F = M/ {^^ J^ ) .M Equation (6) 

The following describes how equations (3) and (4) are obtained. 

Note that the dip vector p represents the apparent dip of a reflector in the 
migrated domain, this reflector passing through the image point x . 

The dip vector p is expressed according to: 

p = (3x3/5x1 , 8x3/8x2, 1 ) Equation (7a) 
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where dxz/dx^ represents the apparent dip along the xi axis and dx^d\2 represents 
the apparent dip along the X2 axis. 

9 and ^ represent the inclination and azimuth respectively of the dipped reflector 
passing through the image point x at which the illumination is evaluated. 

Thus, the following expressions are obtained: 

cos9 = (1+ pi + p2 )"^^^ Equation (7b) 

tan^ = P2 / Pi Equation (7c) 

As we have seen above, starting from the Taylor series development of the 

diffraction travel time (known from the isochronic migration map), the following 

equation is applicable: 

h i^r ) = ix) + (V/, )'' (^, - 3c) H- 1 / 2(x, - xY A^ J, • ( J, - J) + . . . 

Equation (2) 

The reflection point x r is firstly estimated along the reflection plane for which 
the diffraction travel time is stationary (condition for specular reflection). 

In order to consider positions only belonging to the dipped plane, a new vector 
Xp is introduced such that: 



'\ 0 ^ 

0 1 



Equation (7d) 

Note that in this case, as required, the two vectors of the matrix M are actually 
perpendicular to the dip vector p . 

Rewriting equation (2) using equation (7d) gives the following expression: 
^/^p) = '</(^) + (V,/,)' M- +l/2-x; ■M'"-A,,,/, -M^; +... 

Equation (7e) 

namely: Tj{Xp) = Tjix) + P -Xp+l/l-xl Fxl +... Equation(7f) 
In this equation (7f) 
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Consequently, the position for which the diffraction travel time t^ix^) is 
stationary is deduced: 

V^(^p)U^P=0, => x^= -F'*-^ => x^=x+M.x^ 

Equation (7g) 

Equations (7g) and (2) are combined to obtain an expression for the reflection travel time 
in equation (4): 

As already mentioned, the illumination fold is incremented during operation 3 as 
a function of the difference between the diffraction travel time and the reflection travel 

time . 

Even if the reflection travel time t^(xAp);s,r) is different from the diffraction 

travel time (x; J, r ) , due to the limited band nature of seismic data, the migration of the 

reflection that took place at x can constructively contribute to the image of the reflector 
with a dip along x . 

The selection between illuminating migration contributions (constructive 
interference) and non-illuminating migration contributions (destructive interference) is 
made using the increment function i{t^,t/,sj) . 

According to one embodiment of the invention, the increment function is 
expressed according to i{t^ ,t/,s,r) = s(t^ (J; 5", r ) - (x^ (p); ^, r)) where: 

- s(t) is a seismic wavelet which may be monochromatic, with limited band 
or infinite band and it can be chosen to be dependent on time and space, 

s(t) is the derivative with respect to time of the seismic wavelet 
s(t) = d[s(t)]/dt. 
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A derivative with respect to time is usually applied to each seismic trace before 
the Kirchhoff migration (in time or in depth). Thus, particularly in order to make the 
studied quantities uniform, the derivative s(t) of the seismic wavelet with respect to time 
is chosen instead of using the seismic wavelets s(t) itself. 

According to another embodiment of the invention, particularly when the 
derivatives of seismic traces with respect to time are not determined, the increment 
function takes account of the seismic wavelet s(t) instead of the derivative s(t) of the 
wavelet with respect to time; the result is then 
iitj J,;sJ) = s(t^ (x; J, r ) - (J, (p);s, r )) . 

A migration amplitude correction term w(jc;j,r) is conventionally calculated 
during the migration, and in particular takes account of a regular acquisition geometry. 

Therefore, this term w(jc;J,r) may advantageously be used in the method 
according to the invention, and act as "preconditioning'' to the correction made by 
determination of the previously described illumination fold. 

The increment function may then be weighted by this preconditioning term, such 
that the illumination fold l(x,p\s,r) related to a (source J, receiver r) pair is 
incremented by: 

The evaluation of the illumination fold in the migrated domain according to the 
invention is optimum since it is applicable to any acquisition geometry. 

Based on an inverse approach, the illumination fold I(3c,p) estimated by the 
invention agrees with the theoretical results established by Bleistein. This is not the case 
for other methods proposed to solve a similar problem. 

As already mentioned earlier, the evaluation of the illumination fold proposed by 
the invention uses only basic quantities available for migration conventionally done by a 
Kirchhoff summation. 

Advantageously, the illumination fold according to the invention may be 
estimated while seismic traces are being migrated. 
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This is not the case for some other proposed methods (Schneider et al., Laurain 
and Vinje) for which the said evaluation is done post-migration. 

According to another aspect, the invention also proposes a method for correction 
of the amplitudes of seismic data recorded during 3D seismic prospecting in order to 
compensate for the effect of non-uniform illumination of subsoil reflectors, characterized 
in that it comprises steps consisting of: 

- estimating the illumination fold lix.p) by the method according to the 
first aspect of the invention described in detail above, 

- using the inverse l~\x,p) of the said ratio as the weighing factor to be 
applied to each of the said seismic data amplitudes. 

According to yet another aspect, the invention also proposes a method for 
selection of one of a plurality of acquisition geometries as a function of the target of 3D 
seismic prospecting, characterized in that it comprises steps consisting of: 

- for each acquisition geometry considered, determining the illumination 
fold I (J, p) estimated according to the first aspect of the invention given in detail above, 

- selecting the acquisition geometry providing the optimum illumination 
fold as a function of the target. 

Finally, it will be noted that the illumination fold I(3c,p) estimated according to 
the first aspect of the invention described in detail above, may also be used after the 
migration to eliminate acquisition pattems, and to deconvolute the migrated image. 
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